Jamming IV: Distribution of volumes and coordination number in jammed matter: 

mesoscopic ensemble 
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We investigate the distribution of the volume and coordination number associated 
to each particle in a jammed packing of monodisperse hard sphere using the meso- 
scopic ensemble developed in Nature 453, 606 (2008). Theory predicts an exponential 
distribution of the orientational volumes for random close packings and random loose 
packings. A comparison with computer generated packings reveals deviations from the 
theoretical prediction in the volume distribution, which can be better modeled by a 
compressed exponential function. On the other hand, the average of the volumes is well 
reproduced by the theory leading to good predictions of the limiting densities of RCP 
and RLP. We discuss a more exact theory to capture the volume distribution in its 
entire range. The available data suggests a plausible order /disorder transition defining 
random close packings. Finally, we consider an extended ensemble to calculate the 
coordination number distribution which is shown to be of an exponential and inverse 
exponential form for coordinations larger and smaller than the average, respectively, 
in reasonable agreement with the simulated data. 
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I. INTRODUCTION 

Jammed matter refers to a broad class of physical 
many-body systems ranging from granular matter to fric- 
tionless emulsions, and colloids. These systems share the 
property that their constitutive particles can be blocked 
in a configuration far from thermal equilibrium when un- 
dergoing a jamming transition. The statistical mechani- 
cal description of these materials is based on the volume 
fluctuations of the system [l[ taken to be the conservative 
quantity instead of energy, as typically done in thermal 
system. Therefore, the probability distribution of the 
volume occupied by each jammed particle is of partic- 
ular interest and many experimental studies have been 
devoted to study them in detail @, H, 0] . 

Recently, a theory of volume fluctuations has been de- 
veloped at the mesoscopic level providing a relation be- 
tween the volume occupied by each particle and its num- 
ber of contacts [5[. Thus, the distribution of particle vol- 
umes in the system is intimately related to the distribu- 
tion of contacts per particle. In this paper, we calculate 
the distribution of (orientational) volumes occupied by 
each particle in a jammed system of monodisperse hard 
spheres and the distribution of coordination numbers by 
following the theoretical formalism of [H , which is in turn 
based on the Edwards statistical mechanics of jamming 

Q. 

The distributions of volumes and contacts in real pack- 
ings represent ensemble averages in the statistical me- 
chanics sense. Therefore, the distributions depend on 
the state of the packings specified by their compactiv- 
ity through a Boltzmann-like probability in the partition 
function. We show that in two limiting cases of zero 
and infinite compactivity (corresponding to the random 
close packing, RCP, and random loose packing, RLP, re- 
spectively) the distributions can be obtained in analyt- 
ical form. Theory predicts that the distribution of ori- 



entational volumes is exponential with a mean volume 
varying with the average coordination number for RLP 
and constant for RCP. The theoretical predictions are 
compared with computer generated jammed packings of 
equal-size spheres for any friction coefficient. We find 
that the theory well reproduces certain features of the 
numerical distributions, but not all, in the entire range 
of volumes. 

The mean value of the occupied volumes is well re- 
produced by the theory. However, we find important 
deviations between theory and simulations for the higher 
moments of the distribution. For instance, simulations 
show a plateau for small volumes while theory predicts 
an exponential. For intermediate values, the exponen- 
tial shape seems to provide a good fit to the simulated 
data and the predicted dependence of the characteristic 
volume on the coordination number follows partially the 
theoretical prediction. However, a higher scrutiny shows 
deviations in the tail of the distribution which is found to 
decay faster than exponential, having a compressed ex- 
ponential tail. We conclude that the full understanding 
of the distribution requires a more precise theory. We 
discuss how to obtain more exact solutions of the volume 
distribution which can capture the behavior in the entire 
ranges of volumes. 

The study of the distribution of coordination numbers 
is performed by generalizing the ensemble proposed in Q 
to include fluctuations in the number of contacts. The- 
ory predicts an exponential decay for large coordination 
number and an inverse exponential for small coordination 
number. Computer simulations well reproduce the pre- 
dictions. Overall, the present paper serves as a critical as- 
sessment of the theoretical predictions of the mesoscopic 
theory towards the development of an exact formulation 
at the microscopic level that could capture the behavior 
in the entire range of volume and coordination number 
fluctuations of jammed matter. 
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FIG. 1: Schematics of the Voronoi volume and the orienta- 
tional volume associated with particle i. The boundary of the 
Voronoi cell (shown in 2d for simplicity) corresponds to the 
irregular pentagon in black which defines Wj Vor . The limit of 
the Voronoi cell of particle i in the direction s is the mini- 
mum of Tij /2 cos 9ij over all the particles in the packing, as 
indicated. This defines the orientational volume W| which 
is the volume of the sphere of radius rij /2 cos 9ij defined by 
the dash red circle in the figure. The Voronoi volume is the 
integration of the orientational volume over s as in Eq. J2). 
Notice that w s , the orientational free volume associated with 
Wi defined in Eq. Q, ranges from zero (when the orienta- 
tional volume coincides with the volume of the central ball, 
that is when the direction s coincides with a contact point) to 
in principle oo for an isolated ball, although the maximum w s 
in a jammed system is, of course, bounded by the free space 
given by the first or second coordination shell. The Voronoi 
free volume, on the other hand, cannot be zero, by definition. 



II. PDF OF THE ORIENTATIONAL FREE 
VOLUMES 

We start by introducing the different definitions of the 
volume associated with each particle necessary to under- 
stand the problem. The starting point is the volume of 
a Voronoi cell associated to each particle. The Voronoi 
tessellation tiles the entire packing is shown in [f| [f| to 
be a good candidate for the volume function of jammed 
matter. The volume function replaces the Hamiltonian in 
thermal systems, and describes the state of the jammed 
packings in the ensemble average in the partition function 
[l|, H . The Voronoi volume for each particle W™ r gives 
rise to the total volume of the system W = J^iLi W™ r , 
when summed up over all the N particles. In terms of 
the relative coordinates of the particles, we have ob- 
tained in [a, [6[ for monodisperse particles of radius R and 
volume V g the following formula for the Voronoi volume: 



where the integration is performed over the direction s 
forming an angle 9.^ with r^ as in Fig. [TJ and cos 9ij = 
s ■ fij . Taking advantage of this integration we can define 
an orientational Voronoi volume, VV* , for a fixed direction 
s, satisfying: 

wr = jj- s <fwids = (w?) s , (2) 

from which we obtain: 

Wf = V g (^= min ^-)\ (3) 
\2R s-f zj >o s ■ rij J 

W* defines the orientational Voronoi volume which is ob- 
tained without the integration over s. 

The average of the orientational volume over s for a 
single particle, Eq. ([3]), is the Voronoi volume, {W") a — 
W™ 1 and the average of the orientational volume over 
many particles for a fixed s is the same as the average of 
the Voronoi volume over the particles: (W/)» = (W™ r )i, 
in the case of isotropic systems. This last property is 
useful since it allows the use of the orientational volume 
to define the ensemble average of the volume fraction 
without resorting to the use of the Voronoi volume which 
contains the average over s and therefore is more difficult 
to treat from a theoretical point of view. We therefore 
promote the use of the orientational volume function W| 
as the fundamental quantity to characterize the state of 
jammed matter instead of W™ r - It is important to note 
that the probability densities of P(Wf) and P(W™ r ) in 
general differ. For instance, as discussed in Fig. Q] the 
orientational free volume can be for instance zero while 
the Voronoi free volume cannot. 

We define the reduced free orientational volume func- 
tion as 

w = V ■ ( 4 ) 

(we drop the subscript i in w s for simplicity of notation). 

In what follows, we provide a theory for the proba- 
bility distribution function of the orientational free vol- 
ume, P(w s ), which is less complex than the full Voronoi 
volume distribution of Eq. QJ. In [6], this distribu- 
tion is obtained under assumption of uniformity in the 
packing, making the theory valid at a mesoscopic level 
of a few particles diameters. This approximation can 
be seen as defining quasiparticles of free volume w s cap- 
turing the behavior at the mesoscopic distance. Under 
this approximation the inverse cumulative distribution 
P> is obtained (see Eq. (20) in 5] and Eq. (41) in 
from where the probability density can be calculated as, 
P( w s ) = d{1 ~ p >\ then 

V / dw s ' 



Wr oI = \ /(min (^-)Vd S . (1) P(w*) = -exp(-—), (5) 

o J \s- rij >Q is ■ rij I w \ w ) 
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where the average value over the particles, 



w 



(w 



w s P(w s )dw s 



(6) 



was found to be directly related to the so-called geomet- 
rical coordination number z as 



w{z) = -, 
z 



(7) 



where k = 2\/3. 
We note that 



v„ 



- 1 



v„ 



- 1. 



(8) 



Therefore, the orientational volume w s captures the be- 
havior of the average volume function and thus can 
be used as the fundamental variable to define the mi- 
crostates of the system instead of the more complicated 
Voronoi volume. 

The distribution of Eq. © is not the distribution that 
one would obtain in real packings (generated either ex- 
perimentally or numerically) corresponding to the ensem- 
ble average of Eq. ([5]). Therefore, further examination 
is required to derive the distribution of orientational vol- 
umes which can be directly compared with real packings; 
their states determined by the compactivity, X. 

Under the volume ensemble point of view [l[ , the ob- 
servables in real packings are ensemble average over the 
Boltzmann distribution function determined by the aver- 
age volume w. The probability distribution of volumes in 
the canonical volume ensemble for a single quasiparticle 
of orientational volume w is then: 



P(w s \X, Z) = ^-J P(w s ) exp (-|) g(w)dw, (9) 

where we have introduced the density of states, g(w), and 
the partition function is 



-Ziso = / exp(—w/X)g(w)dw. 



(10) 



Under the considerations of |5j, the volume function 
depends on z, Eq. ([7]), and we have an exponential 
density of states, g(z) = (h z ) z - 2d = e - ^ -2 ^/**, with 
the constant h z <C 1, representing the typical separation 
of the configurations in the phase space (analogous to 
the Planck constant in quantum mechanics). We have 
z* = —X/]nh z . The limit of integrations in ((§]) are given 
by the isostatic condition 0, d, HI over the mechanical 
coordination number, Z (then, the subscript iso in the 
partition function, Eq. (|10p h counting the contacts with 
nonzero forces and is different from the geometrical co- 
ordination number, z, counting all contacts, even those 
with zero forces. We have Z < z < 6 [5]. The mechan- 
ical coordination, Z(fi), depends on friction and varies 
between Z(0) = 2d = 6 and Z(oo) = d + 14. 
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FIG. 2: Computer generated packings arranged in the phase 
diagram of [Ej. The packings are used to calculate the dis- 
tribution of orientational volumes and coordination number. 
We concentrate our study on packings generated with differ- 
ent friction as indicated in the figure following the RLP line 
from the J point to the L point and the RCP line from the J 
point to the C point. 



The distribution of orientational volumes for a packing 
at a given compactivity is then: 



1 M*= z ) i / , r 
P(w a \X,Z) = — I -exp(- — 



--iso Jw(z=6) 



exp -- 



2V3 



(11) 



In h, dw. 



This equation cannot be solved analytically for a gen- 
eral (X,Z). However, analytical forms can be obtained 
in the limiting cases of X = (the ground state) and 
X — ► oo: the RCP and RLP lines in the terminology of 
respectively. The advantage of studying these dis- 
tributions is that they can be checked with simulations 
or experiments without the use of the compactivity as a 
fitting parameter. 

From Eq. (fTTj) . we find along the RCP line, 
Prcp(w s \Z) = P{w s \X = 0,Z): 

Prcp(w s \Z) = V3cx P ^-w s V3), 2 6(4,6], (12) 
and for the RLP line, P RL p(w s |Z) = P( W S \X -> oo, Z): 

PRLPK|Z) = -^=exp(-^|), Ze[4,6]. (13) 



We note that both limiting distributions coincide at Z = 
6, the frictionless J-point. 

In what follows, we test the above predictions with 
computer simulations. We generate packings at the jam- 
ming transition using the split algorithm explained in 
The packings consist of 10,000 spherical equal-size 
soft particles interacting via Hertz normal forces, Mindlin 
tangential forces and the Coulomb condition with fric- 
tion coefficient [i. The mechanical coordination number 
Z versus the volume fraction <j> of the generated packings 
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FIG. 3: PDF of w 3 . Exponential fit in a semilog plot, (a) 
We plot the results for the distribution of the packings along 
the RCP line with different values of friction, as indicated in 
the figure (see Fig [2J ■ The red dashed line is a fit with the 
theoretical prediction Prcp(w s ), Eq. (112[) . but with the in- 
verse characteristic volume or slope of 2\/3 instead of %/3~ as 
predicted by the theory, (b) Same as (a) but for the packings 
along the RLP line in Fig. [2] prepared with different /i. The 
red and blue dashed lines are fits to the theoretical predic- 
tion -Prlp(w s ), Eq. (|13|) . but with the inverse characteristic 
volume or slope replaced from Zj{2\pS) — > Z/y/3, twice as 
predicted by the theory like in (a). The red line is for /i = 
ad has slope 2y/3. The blue line is for fi — > oo and has slope 

4/V5. 



are plotted in Fig. [2] in the framework of the phase dia- 
gram of [5] . We change friction from /j, = to fx — *■ oo to 
generate the packings along the RCP line as indicated in 
the figure with the corresponding change in the mechan- 
ical coordination number from Z(0) = 6 to Z(oo) = 4. 
The RCP line is also generated by changing friction but 
the volume fraction remains constant, as seen in the fig- 
ure, while the mechanical coordination varies from 6 to 
4. 

We focus on the calculation of the probability distribu- 
tion function of w s for the packings along the RCP-line 
to test Eq. ([12]) and along the RLP-line to test Eq. (fl~3|) . 
Figure [3] shows the results. Along the RCP-line, shown in 
Fig. [3Ja,, we find all distributions are the same, indepen- 
dent of friction and Z, as suggested by Eq. (fl~2"j) . On the 
other hand, the distribution along the RLP line, shown 
in Fig. [3Jd, depends on friction and therefore on Z(/i) as 
suggested by theory. The exponential dependence seems 
to be captured upon a first inspection of the data done in 
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FIG. 4: PDF of w s . Compressed exponential fit in a dou- 
ble log plot. The results are the same as in Fig. [3] but now 
reploted in a double log plot to obtain the compressed expo- 
nential fitting exponent /3 which is extracted from the sloe of 
such a plot, as explained in the text, (a) PDF for the packings 
along RLP-line. (b) PDF for the packings along the RLP-line. 



a semi- log graph of Figs. [3^, and b (arguably better for 
the RCP case), at least for intermediate values and the 
tail of the distribution. In the case of the RLP line (Fig. 
|3b) , the slope of the semi- log plot of the exponential fit 
has the same dependence on Z as predicted by theory 
(that is the slope in the semi-log plot increases linearly 
with Z) but with twice the constant value as predicted 
by Eq. (flU)) . We find that the exponential fit leads to a 
tail with characteristic volume = Z/-\/3, twice the value 
predicted by theory z?/(2v / 3), Eq. (fl3|) . However, the 
linear trend with Z is observed in the data. 

On the other hand, the average of the distributions 
agrees very well with simulations (see below). But when 
we fix the mean according to theory the exponential tail 
is inaccurate. Theory either provides the exponential fit 
with the incorrect average or provides the correct average 
with the deviations from the exponential fit. 

The same situation is observed in Fig. [3^, for the RCP 
line. We force the fitting to be exponential, and then 
the average value has to be modified from V3 — > 2^3- 
The reason for this discrepancy is that the theory does 
not capture the distribution in the full range of volumes. 
Figure [3] clearly show a plateau at small values of w s 
deviating from the exponential behavior predicted by the 
theory. 



5 



Furthermore, a more strict scrutiny of the data seems 
to indicate that the exponential fit may not be sufficiently 
accurate as shown in Figs. [3^ and b. While tempting to 
conclude that the exponential is a good fit to the data 
(for instance, the fitting in Fig. [3^ looks convincing), 
further scrutiny reveals important deviations in the tail. 
To visualize the deviations, one should take the plot and 
look at it, not frontally, but from the side. It is evident 
that there is a slight curvature in the distributions devi- 
ating from the linearity in the semi-log plots. Indeed, the 
distributions decay slightly faster than the pure exponen- 
tial decay predicted by theory. The largest evidence of 
this deviation is perhaps in the tail of the fj, — ► oo data 
in the RLP, Fig. H>. 

A double log analysis of the data shown in Figs. [4^ 
and HJa reveals that a compressed exponential behaviour 
might better capture the tail of the distributions above 
the average value: 



P(w s \Z) ~ Aexp [- (—) " 
L \w c J 



w > w, 



(14) 



where f3 w « 1.5 is the compressed exponential exponent 
(J3 W = 1 would be a pure exponential) valid for all the 
RCP and RLP packings according to Fig. [4^ and|4jD, w c 
is a characteristic volume independent of Z in the RCP 
packings and depends on Z for the RLP packings, and A 
is a constant. 

We want to stress the difficulties associated with a fit 
to a compressed (or stretched) exponential function like 
Eq. ([T4"|) . A double log plot gives: 



In 



-\n(P{w s \Z)/A)) = f3 w \n{w s )- (3 w ln(w c ), (15) 



providing a linear fit with slope (3 W . Beyond the inherent 
subtleties associated to taking a double log of a function 
in a such a short range, a further complication arises 
because such a linear fit depends on the value of the 
constant A. Since Eq. (fl"4"|) covers only the tail of the 
distribution, A cannot be obtained from normalization, 
remaining as a fitting parameter. The result is a fitting 
of f3 w dependent slightly on A, providing an extra level 
of difficulty. 

Although the compressed exponential seems to fit the 
data better than the pure exponential, before more the- 
ory or numerical/experimental evidence in the limit N — > 
oo become available, we are inclined to conclude that the 
problem is not closed. We note that a similar dichotomy 
between exponential and compressed/stretched exponen- 
tial behavior has plagued the study of the distribution of 
forces in jammed matter since the first studies on the 
subject. Given the inherent difficulties in any numerical 
estimation, the dispute will have to eventually be settled 
when more exact theories become available. 

Beyond the distributions of volumes, the theory repro- 
duces very well the average value of the volumes, Eq. [[7]). 
Based on the properties of the averages expressed above, 
there is no need to calculate the full Voronoi volume to 
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FIG. 5: Relation between local volume per particle and its 
coordination number for (a) the volume function and (b) the 
volume fraction. We use the packings along the RLP line with 
the friction as indicated and we calculate the volume function 
binning the data by the coordination number of the particles. 
The red dots correspond to the average over all the particles 
in the packing of the volume function which fits the theory 
very well. The red and blue dashed lines are logarithmic- 
like fittings to the data and we add the value at z = 12 for 
the FCC and the dodecahedron. The data from the packings 
goes only up to z = 11, indicating the absence of ordered 
structures. 



obtain the average volume fraction, since the average of 
the orientational w s suffices. For instance in the friction- 
less packing we find (w s ) — 0.561, which gives a volume 
fraction <f> = 1/(1 + (w s )) = 0.641, in agreement with the 
direct measurement of the volume fraction of the packing, 
0.64. 

A full comparison between theory and simulations is 
given in Fig. [5l where we study the dependence between 
average volume and coordination number. For each pack- 
ing along the RLP line wc calculate the average orienta- 
tional volume focusing on the particles with a given z. 
We also calculate the average over all the particles for 
a given packing, plotted as the red dots in Fig [5] In 
practice, we do not measure the geometrical coordina- 
tion number z but the mechanical coordination number 
Z . However, we know that for the packings along the 
RLP line z w Z j| (this is because h z — > 0). Further- 
more, the RLP line corresponds to X — ► oo and therefore 
the prediction of the average volume function, Eq. ([7]) 
can be tested directly with these fully random numerical 
packings, extending this result, valid for a quasiparticle, 
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to the entire packing. Thus, the packings along the RLP 
line reveal the approximate behaviour of quasiparticles of 
fixed coordination z. We notice that there could be still 
some subtleties when comparing Eq. ([7]) , valid for quasi- 
particles, to the results in packings. Using </> _1 = w + 1, 
we plot the volume fraction in Fig. [SJd. 

Figures [5^ shows that the mean of the distribution of 
volumes, w, is well captured by the theory of Eq. ([7]) 
(see the black dashed line in comparison with the red 
dots in Fig. [5]). This is why the theory provides very 
good fittings to the values of RCP and RLP in [5] . The 
agreement can be seen as well in the volume fraction in 
Fig. 03, and exists despite the fact that the full distri- 
bution presents the deviations discussed above. 

Figure [5] presents further interesting results. For a 
given packing along the RLP-line specified by a fixed 
friction, there are a variety of particles with varying co- 
ordination z, following a well-defined functional relation 
between the volume occupied by the particle and its coor- 
dination number (see for instance the red and blue dashed 
lines in Fig. [5^ corresponding to fittings for the cases 
(X = and fj, — > oo, respectively). While this plot does 
not tell us how many particles there are for a given coor- 
dination number (see next section) we see that for each 
packing, there exists a variety of local volume functions 
with z ranging from z = (since there are some rattlers) 
up to z = 11 (but not 12, interestingly, see below). 

The assumptions for the limits of integration in the 
partition function in [j| or the ensemble average Eq. (jlip , 
Z < z < 6, seem to be violated here. However, we have 
to remember that the theory is mesoscopic and further 
coarse-graining is needed for these bounds to be more ac- 
curate. Regardless, even though the range in z extends 
further than the bounds, Fig. [5] corresponds to the av- 
erage for a fixed z but does not tell how many states 
there are for every z. When these details are properly 
taken into account, the bounds are approximately sat- 
isfied, although fluctuations persist, bringing us to the 
next Section of this paper. 

Focusing around the z = 12 point in the figures, we 
observe that an extrapolation of the fitting to the curve 
of frictionless packing seems to converge to the green dot 
in Figs. [5] at z = 12 which correspond to the free volume 
function of FCC, w FCC = 0.35135 (Fig. [5k,) and the FCC 
volume fraction FCC = tt/\/T8 « 0.7402 (Fig. Eb). An 
extrapolation of the fitting to the infinite friction data 
seems to pass through the volume fraction of the dodec- 
ahedron as indicated by the blue point in Fig. [5J which 
has also z = 12 but slightly larger volume fraction that 
FCC [the dodecahedron can't tile the space without leav- 
ing holes, so the best global packing is still the FCC]. We 
observe that there are no particles in the packings with 
z = 12. Indeed the green and blue points at z = 12 in 
Figs. Ek an d b were added by hand and the real curves 
stop shortly at z — 11. The absence of z = 12 states 
indicate the randomized state of the systems. 

More importantly, we see that if we extend the theo- 
retical result of Eq. ([7]) to the ordered region, examining 



z from z = 6 to 2 = 12, the theory does not fit the FCC 
value. Instead, we obtain w(z = 12) = 2^3/12 = 0.2886, 
below the FCC or dodecahedron value. In principle, 
this result is expected since the theory assumes random 
isotropic states while the FCC is an ordered anisotropic 
packing. However, the absence of a good fitting of the 
disordered branch through the FCC, together with the 
fact that the theory fits so well the disordered states, 
raises the interesting question of the existence of a phase 
transition between the RCP limit at z = 6 and the FCC 
at z = 12. Since packings cannot equilibrate above z = 6 
without the formation of crystalline regions, we expect 
an ordered branch from the FCC point towards the RCP 
point. It seems plausible that there could be a disconti- 
nuity from the disordered branch to the ordered branch, 
the existence of which could determine whether there ex- 
ist a disorder/order phase transition characterizing the 
RCP; a question of interest in the field. 

To summarize this part of the study, while theory 
predicts an exponential behavior and approximates well 
some features of the distribution such as the average 
value, simulations indicate that a compressed exponential 
fitting could be also possible. We therefore require more 
refined theories to account for the full behavior of the 
volume distribution. The present approach suggests that 
study of the ensemble of quasiparticles of fixed coordina- 
tion number could provide clues to the behavior of the 
entire system through the Edwards ensemble approach. 
While new theoretical concepts are required, current at- 
tempts indicate that it might be possible to develop a 
theory of the volume distribution that is exact to, at 
least, a given coordination shell of particles. It appears 
that P(w s ) for a fixed z-ensemble might be solved exactly 
by a brute force approach, since the range of a Voronoi 
cell is finite. Although, it may contain a large number of 
variables, the computer should handle such a computa- 
tion. Such an analysis parallels the Hales proofs of the 
Kepler's conjecture [10| . 

III. PDF OF THE COORDINATION NUMBER 

Next, we analyze the distribution of coordination num- 
ber by generalizing the theory of [5(]to include fluctua- 
tions in z. While we have assumed |5| that every single 
quasiparticle satisfies the specified bounds: Z < z < 6, 
below we relax this constraint to extend the bounds to 
the geometrical coordination of the entire system by con- 
sidering: 

N 

Nz min <y^ Zj< Nz max , (16) 

i=l 

where in the following we set z m j n = Z and z max = 6. 
This new condition implies that is not possible to con- 
sider the single quasiparticle partition function, Eq. 
used in || and that the full TV-particle partition func- 
tion has to be considered. Then, the ensemble average of 
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FIG. 6: P(z): comparison between theory and simulations in 
(a) a lin lin plot and (b) a semi-log plot to appreciate better 
the tail of the distribution. We plot the systems at the J-point 
(jj, = 0) and L-point (/i = oo) in Fig. [2] 



the probability distribution function of the geometrical 
coordination number, P(z), is 



/v 



1 

"z 



(JV2min-z)<Ei=l z z <(Nz max -z) 



(17) 



N-l 



i=i 

where the partition function is 



Z = 



N 



TX e -(z</z*+Mz<) dZi . 

'jVZmin^i:™ 1 2<<JV* m „ fj{ 

(18) 

We have z* = —l/\nh z , the inverse compactivity (3 — 
1/X, and the coupling constant 



(19) 



We obtain: 



P(z) 



Erf( v^ES^zM) + Erf( ^B^J ) 



Erf( VJV ^:^ ) + Erf( ^^W ) 

~{z/z*+p K /z) 



(20) 



with the following constants: 



z)/(N - 1) « z max + (z max - z)/iV 

/V. 
(21) 



z min = (iVz min - z)/(7V - 1) « z min + (z min - z)/N 



The constants /I and a in Eq. (|2H)) are significant be- 
cause they represent the mean and standard deviation 
of the a Gaussian expansion in a saddle-point approxi- 
mation of the inverse Fourier transform of the partition 
function allowing the calculation of the free- volume den- 
sity. They are: 



. 1/2 K X {2B^) 



and 



\K (2B 1 / 2 ) K (2B 1 / 2 ) 2 



(22) 



(23) 



where K n {a) is the modified Bessel function of the second 
kind: 

/>oo 

/ x n e-% (x+1/x) dx = 2K n {a). (24) 
Jo 

Furthermore, we approximate: 

(n - i)(z; iax - pf w(jv - i)(( 2max - pf+ 

+ 2(z max - p)(z max - z)/N) 
~N{z max - fl) 2 + 

+ (Zmax - /i)(Zmax +jX~2z), 



(25) 



and 



(N-l)(z' m 



p) 2 ^(N-l)((z mia -p) 2 + 

+ 2(z mhl - p){z min - z)/N) 

~N(z min - pf + 

+ (fi- z min )(2z - z min - p). 

(26) 

Next, we consider two cases: When z max — p < p— z min , 
then p > (z max + z m i n )/2, and we obtain: 



P(z) oc exp 



* 1 ~ Z 2 

z* \ u 2 



(3k 
z 



(27) 



8 



Otherwise, we obtain 
P(z) oc exp 



~~ V 2 j ~~ 



(28) 



We consider the approximation: 



ct 2 £ + L(B) - L(B) 2 

B \ 1/2 e 

-g- , - j . 



(29) 



where e is the ensemble average of the volume function 
(analogous to the energy density): 



n e dB ' 
which is related to the free volume density /: 

^ TV^oo AT 

We consider Z max = z max /z*, Z min = z min /z* and 



r/m- 7?1 / 2 ^l(2BV 2 ) 
L(S) = i? K (2Bi/2)' 



(30) 



(31) 



(32) 



and we find that B^ 1 / 2 ps ez*/n 1 (3 — Bz*/k = k/ (z*e 2 ), 
and the maximum volume density is e max w L(B)/f3 — 
k/z*B^ 1/2 = K/(z*Z min ) = n/zmin. Similarly e m i„ sa 

^/ 2max- 

With Eq. (JUJ), the distributions Eqs. {57]) and 

become: 



and therefore D/C — z 2 lax — 36, z max = 2c? = 6. While 
the constant z* (which determines the density of states) 
is difficult to know a priori, the ratio D/C can be eas- 
ily compared with the simulations since it relates to the 
mechanical coordination number only and not z* . 
For the RLP L-point we find: 



Prlp(z) tx exp 





z 2 M 

mm \ 







(37) 



= d + 1 = 4, 



P(z) oc exp 
and 



and therefore D/C = z 2 lin = 16, z n 
where again the ratio is given in terms of the average 
coordination number, this time in terms of the lower limit 
as expected for the RLP state. 

We test these forms with the computer generated pack- 
ings at the J-point and L-point. Figure [6] shows the re- 
sult. The lin-lin plot of Fig. [S^ shows that the distri- 
bution near the average value is well approximated by 
the theory for both points. To investigate the tails of the 
distributions, Fig. [6Jd plots a semi-log curve. While some 
deviations are observed, the fit is still reasonable except 
for the larger coordination number which are very rare, 
3 x 10~ 2 less probable than the most probable value. 
The low z behavior corresponds to the dependence com- 
ing from the Boltzmann distribution of volume while the 
large z dependence depends on the shape of the density of 
states. Thus, the deviations at large z can be attributed 
to deviations from the assumed exponential form. For 
coordinations from 2 to the theory deviates from the 
data, but these contacts corresponds to rattlers. 

The coefficients of the best fits agree with the theory. 
We find for the J-poi nt: C — 1.55 and D = 52.2, with 
D/C = 33.68, and ^/DjC = 5.8, 3% off the value pre- 
dicted by theory, yj D/C = 6. For the RLP L-point we 



L( z (2-e/e min ) + ~)], (33) find C = 1.25 and D = 18.10 with \J~D~JC = 3.81, close 
V 2 e ) J to the prediction \JD/C = 4. 



P(z) oc exp 



z (2 - e/e max ) + 



z* \ z e 

Therefore, the probability of z has the shape: 



P(z) oc exp 



(34) 



(35) 



where C and D are constants independent of z and de- 
termined by the ensemble average volume function, e, 
and the compactivity, X . At the frictionless J-point and 
the infinitely frictional L-point (see Fig. [2|) these con- 
stants take simple forms and the distributions reduce to 
the following. 

For the RCP J-point we find: 



-Prcp(-z) oc exp 



z + 



(36) 



IV. CONCLUSIONS 

We have presented the predictions of the mesoscopic 
theory presented in Q concerning the probability dis- 
tribution of the orientational volumes in jammed matter. 
The theory captures very well the average volume and in- 
deed gives rise to good predictions of the RCP and RLP 
volume fraction as shown in Q . However, when compar- 
ing the full distribution we find important deviations. For 
instance, computer simulations are able to detect slight 
deviations from the pure exponential behaviour predicted 
by theory. This deviation could be better approximated 
by a compressed exponential behavior, although a more 
conclusive fitting necessitates a more precise theory or 
simulations in the large scale limit. The plateau observed 
in the volume distribution is not captured by the the- 
ory either. The present approach indicates that further 
progress could come from a systematic analysis of the 
quasiparticles in a fixed z-ensemble. 
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On the other hand, the behavior of the probability dis- 
tribution of coordination number is captured well by the 
theory. Here, the mesoscopic theory of considering 
restricted bounds for each quasiparticle, is extended to 
allow for coordination number fluctuations, by imposing 
the bounds to the entire system. Such a problem can 
be solved under several approximations and predicts a 
mixed exponential forms that are well reproduced by the 
simulations at the J and the L-point. Such an extended 
z-ensemble can be also applied to obtain predictions of 
the phase diagram, a study that provides interesting re- 
sults, such as a possible phase transition at RLP, but will 
be reported in a follow up paper. 



Why does the theory capture the contact distribution 
better than the volume distribution? The distribution of 
contacts is an ensemble average based on the states char- 
acterized by w(z). Since this function is quite accurate, 
the contact distribution follows. On the other hand, the 
volume distribution is based on the uniform approxima- 
tions done on Q, and therefore is not quite exact. The 
full distribution of volumes is where more theoretical de- 
velopments are required. 
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